1 REM ================================= 2 REM 3 REM DIOPHANTINE SOLVER 4 REM RUPERT REPORT #62 5 REM --------------------------------- 6 REM SOLVE EQUATIONS OF THE FORM 7 REM AX + BY = C FOR X AND Y 8 REM WHERE ALL TERMS ARE INTEGERS 9 REM ================================= 10 GOSUB 100 :REM INITIALIZE 20 GOSUB 300 :REM GET A,B,C 30 GOSUB 1000 :REM MAKE A & B POSITIVE 40 GOSUB 2000 :REM CONVERT A/B TO CONTINUED FRACTION AND GET GCF(A,B) 50 GOSUB 3000 :REM TABLE OF CONVERGENTS 60 GOSUB 4000 :REM FIND T FOR X,Y > 0 70 GOSUB 5000 :REM FIND X,Y FOR GIVEN T 90 END 99 REM -------------------------------- 100 REM -- INITIALIZE & INPUT VARIABLES 110 DIM P(20),Q(20),F(20) 120 DEF FNDIV(B)=INT(A/B) :REM INTEGER DIVISION 130 DEF FNMOD(B)=INT(A-B*FNDIV(B)) :REM MOD FUNCTION 140 DEF FNII(Z)=(Z=INT(Z)) :REM TRUE IF Z IS INTEGER 150 DEF FNMI(X)=SGN(X)*INT(ABS(X)) :REM MAKE INTEGER 160 DEF FNX(T)=C*Q0-SNB*B*T :REM X(T) 170 DEF FNY(T)=A*T-SNB*C*P0 :REM Y(T) 175 PRINT "[147]" 180 PRINT: PRINT "ENTER NON-0 INTEGER COEFFICIENTS A,B,C" 190 PRINT "FOR EQUATION AX + BY = C:" 200 INPUT A,B,C 210 IF NOT (FNII(A) AND FNII(B) AND FNII(C)) THEN PRINT"INTEGERS ONLY": GOTO 180 220 IF A=0 OR B=0 OR C=0 THEN PRINT"NON-ZERO INTEGERS ONLY": GOTO 180 230 PRINT "[147]EQUATION IS ";A;"* X + ";B;"* Y = ";C 240 PRINT 250 A0=A: B0=B: C0=C: REM SAVE ORIGINALS 260 RETURN 300 REM -- CONVERT TO POSITIVE INTEGERS 310 SNB=1 :REM SIGN OF B 320 IF A<0 THEN A=-A: B=-B: C=-C :REM MAKE A>0 330 IF B<0 THEN B=-B: SNB=-1 :REM MAKE B>0 340 A1=A: B1=B: C1=C :REM SAVE WORKING VALUES 350 RETURN 360 REM -------------------------------- 1000 REM CONVERT A/B TO CONTINUED FRACTION F(N) AND FIND GCF 1010 REM A/B = F(N) REMAINDER M 1020 N=1 :REM NUMBER OF TERMS 1030 F(N)=FNDIV(B) :REM INTEGER DIVISION 1040 M=FNMOD(B) :REM REMAINDER 1050 IF M=0 THEN 1080 :REM DONE 1060 A=B: B=M: N=N+1 1070 GOTO 1030 1080 IF FNII(N/2) THEN 1110 :REM N/2 IS INTEGER::N IS EVEN 1090 F(N)=F(N)-1 :REM CONVERT TO EVEN NUMBER OF TERMS 1100 F(N+1)=1: N=N+1 1110 GCF=B :REM FINAL DIVISOR IS GCF OF A,B 1120 A=A1: B=B1: C=C1 :REM RESTORE WORKING VALUES 1130 REM -----NEXT LINE PRINTS CONTINUED FRACTION FOR A/B 1140 PRINT A; "/"; B; "= [";: FOR K=1 TO N: PRINT F(K);: NEXT K: PRINT "]" 1150 PRINT " GCF ("; A; ","; B; ") ="; GCF: PRINT 1160 RETURN 1170 REM ------------------------------- 2000 REM --REDUCE EQN BY GCF IF POSSIBLE 2010 IF GCF=1 THEN 2050 2020 IF FNII(C/GCF) THEN GOTO 2040 2030 PRINT "NO SOLUTIONS. C IS NOT DIVISIBLE BY GCF OF A AND B": END 2040 A=A1/GCF: B=B1/GCF: C=C1/GCF 2050 RETURN 2060 REM ------------------------------- 3000 REM -- TABLE OF CONVERGENTS OF A/B 3010 P(0)=0: P(1)=1 3020 Q(0)=1: Q(1)=0 3030 REM PRINT "CONVERGENTS P,Q:" 3040 FOR K=2 TO N+1 3050 P(K)=F(K-1)*P(K-1)+P(K-2) 3060 Q(K)=F(K-1)*Q(K-1)+Q(K-2) 3070 REM PRINT P(K),Q(K) 3080 NEXT K 3090 P0=P(N): Q0=Q(N) 3100 RETURN 3110 REM ------------------------------ 4000 REM -- SHOW RESULTS 4010 PRINT: PRINT "X AND Y AS FUNCTIONS OF T :" 4020 PRINT " (T = 0, +/-1, +/-2, ...)": PRINT 4030 SN$=") - (": IF SNB=-1 THEN SN$=") + (" 4040 PRINT " X = ("; C*Q0; SN$; B; "* T )" 4050 PRINT " Y = ("; A; "* T "; SN$; C*P0; ")" 4060 PRINT: PRINT "VALUES OF T FOR X & Y > 0:" 4070 TX=SNB*C*Q0/B 4080 TX$=" < ": IF SNB=-1 THEN TX$=" > " 4090 TY=SNB*C*P0/A: TY$=" > " 4100 PRINT " FOR X>0: T"; TX$; TX 4110 PRINT " FOR Y>0: T"; TY$; TY 4120 IF SNB=-1 THEN 4200 4130 REM FIND T FOR TTY 4140 IF TX < TY THEN PRINT "NO INTEGER SOLUTIONS FOR X>0 AND Y>0": GOTO 4230 4150 TA=INT(TX): IF TA=TX THEN TA=TA-1 4160 TB=INT(TY)+1 4170 IF TATX AND T>TY 4210 MN=INT(TX)+1: IF TY>TX THEN MN=INT(TY)+1 4220 PRINT "T IS ANY INTEGER GREATER THAN OR EQUAL TO "; MN 4230 PRINT 4240 RETURN 4250 REM ------------------------------ 5000 REM --FIND X,Y FOR GIVEN T 5010 INPUT "WHAT VALUE OF T";T 5020 PRINT " X IS ";FNX(T) 5030 PRINT " Y IS ";FNY(T) 5040 PRINT: INPUT"MORE [Y,N]Y[157][157][157]"; A$ 5050 IF A$="N" THEN RETURN 5060 PRINT "[145] ": PRINT"[145]"; 5070 GOTO 5010